“Calhoun 


Institutional Archive of the Naval Postgraduate School 





Calhoun: The NPS Institutional Archive 
DSpace Repository 


Theses and Dissertations 1. Thesis and Dissertation Collection, all items 


2003-06 


Determination of bulk mechanical properties 
of nanostructures from molecular dynamic simulation 


Duff, Richard A. 


Monterey, California. Naval Postgraduate School 
http://hdl.handle.net/10945/994 


This publication is a work of the U.S. Government as defined in Title 17, United 
States Code, Section 101. Copyright protection is not available for this work in the 
United States. 


Downloaded from NPS Archive: Calhoun 


Calhoun is the Naval Postgraduate School's public access digital repository for 


\§ D U DL EY research materials and institutional publications created by the NPS community. 
«iit Calhoun is named for Professor of Mathematics Guy K. Calhoun, NPS's first 


NNN KNOX appointed -- and published -- scholarly author. 


LIBRARY Dudley Knox Library / Naval Postgraduate School 
411 Dyer Road / 1 University Circle 


http://www.nps.edu/library Monterey, California USA 93943 


NAVAL POSTGRADUATE SCHOOL 
Monterey, California 





THESIS 


DETERMINATION OF BULK MECHANICAL PROPERTIES OF 
NANOSTRUCTURES FROM MOLECULAR DYNAMIC SIMULATION 


by 


Richard A. Duff 





June 2003 


Thesis Co-Advisors: Young W. Kwon 
James H. Luscombe 





Approved for public release; distribution is unlimited 





TH 


IS PAGE 








NTENT 





ONALLY LEFT BLANK 





REPORT DOCUMENTATION PAGE 


Public reporting burden for this collection of information is estimated to average 1 hour per 
response, including the time for reviewing instruction, searching existing data sources, 
gathering and maintaining the data needed, and completing and reviewing the collection of 
information. Send comments regarding this burden estimate or any other aspect of this collection 
of information, including suggestions for reducing this burden, to Washington headquarters 
Services, Directorate for Information Operations and Reports, 1215 Jefferson Davis Highway, Suite 
1204, Arlington, VA 22202-4302, and to the Office of Management and Budget, Paperwork Reduction 
Project (0704-0188) Washington DC 20503. 
1. AGENCY USE ONLY ‘ REPORT DATE 3. REPORT TYPE AND DATES 
(Leave blank) June 2003 COVERED 

Master’s Thesis 
4. TITLE AND SUBTITLE Determination of Bulk Mechanical . FUNDING NUMBERS 
Properties of Nano Structures Using Molecular Dynamic 
Simulation 
6. AUTHOR (S) Richard A. Duff 
7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS (ES) 8. PERFORMING ORGANIZATION 
Naval Postgraduate School REPORT NUMBER 
Monterey, CA 93943-5000 


9. SPONSORING / MONITORING AGENCY NAME(S) AND 10. SPONSORING/MONITORING 
ADDRESS (ES) AGENCY REPORT NUMBER 


11. SUPPLEMENTARY NOTES The views expressed in this thesis are those of the author 
and do not reflect the official policy or position of the U.S. Department of Defense 


or the U.S. Government. 
12a. DISTRIBUTION / AVAILABILITY STATEMENT 12b. DISTRIBUTION 


Approved for public release; distribution is unlimited CODE 





13. ABSTRACT (maximum 200 words) 


Determining bulk mechanical properties from microscopic forces has become 
important in the light of utilizing nano-scale systems. The molecular dynamics model 
was used to determine the modulus of elasticity and shear modulus of pure metallic 
micro lattice structures. Preliminary results indicate that the modulii of elasticity 
is determined to within 15% accuracy for 5 different metals of 500-atom structures 
when compared to the experiment values of bulk materials. Furthermore the elastic 
modulus for copper structures was computed with different temperatures, different 
magnitudes of stresses and various kinds of dislocations. From the preliminary 
results, it is concluded that the model accurately determines the mechanical 
properties of the nano-scale systems. 














14. SUBJECT TERMS Molecular Dynamics, Paradyn, Quantum Mechanics, 
Modulus of Elasticity, Shear Modulus, Modulus of Rigidity 





16. PRICE CODE 
17. SECURITY 18. SECURITY 19. SECURITY 20. LIMITATION 
CLASSIFICATION CLASSIFICATION OF THIS CLASSIFICATION OF OF ABSTRACT 
OF REPORT PAGE ABSTRACT 
Unclassified Unclassified Unclassified UL 


NSN 7540-01-280-5500 Standard Form 298 (Rev. 2-89) 
Prescribed by ANSI Std. 239-18 








TH 


IS PAGE 








NTENT 





ONALLY LEFT BLANK 





el: 


Approved for public release; distribution is unlimited 


DETERMINATION OF BULK MECHANICAL PROPERTIES OF NANO 
STRUCTURES FROM MOLECULAR DYNAMIC SIMULATION 


Richard A. Duff 
Lieutenant Commander, Canadian Navy 
B. Eng., Technical University of Nova Scotia, 1985 








Submitted in partial fulfillment of the 
requirements for the degree of 





MASTER OF SCIENCE IN PHYSICS 


from the 


NAVAL POSTGRADUATE SCHOOL 





June 2003 
Author: Richard A. Duff 
Approved by: Young W. Kwon 


Thesis Co-Advisor 


James H. Luscombe 
Thesis Co-Advisor 











William B. Maier 
Chairman, Department of Physics 








aba 





TH 


IS PAGE 








NTENT 





ONALLY LEFT BLANK 





iv 


ABSTRACT 





Determining mechanical properties from microscopic 





forces has become important in the light of utilizing nano- 


scale systems. The molecular dynamics model 








was used to 


determine the modulus of elasticity and shear modulus of 














pure metallic micro lattice structures. Preliminary 








results indicate that the modulii of elasticity is 








determined to within 15% accuracy for 5 different metals of 
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XV1 


Es INTRODUCTION 


A. OBJECTIVES 





As technology advances, there ar ver-increasing 
demands to place more and more technological advances into 


the same volume. This volume limitation is especially true 





for systems that operate in unfriendly environments, for 
example space. A shuttle that operates in space is limited 
by its carrying capacity. Equipment that is necessary for 


monitoring shuttle systems currently is typically hard 





mounted on surfaces and takes up precious bulkhead area. 








As technology develops, there is an ever-increasing 


pressure to miniaturize such instruments, and perhaps embed 





them into the structure. A key issue is whether the 


structure can retain its material properties and not be 





weakened by th mbedded device. Carbon nanotube materials 





have shown excellent properties such as elasto-mechanical 


and both semi-conducting or metallic, depending on their 








configuration. [1] Properties which give carbon nanotubes 
excellent prospects for embedded nano-scale sensors. 


Carbon nanotubes have also shown that they can be formed 





into various components from the nano-transistor to nano- 





mechanisms, such as pumps and containment devices[2]. 





Embedding carbon nanotube structures and sensors would be 





beneficial in reducing the real estate costs normally 











associated with carrying externally mounted devices. This 
thesis is a provisional investigation into the mechanical 
properties of nano-scale systems embedded in bulk 





materials. We examine the link between microscopic and 
macroscopic behaviors of such materials. Specifically, 


this thesis will study the modulus of elasticity(E) and 





1 


shear modulus(G) that can be obtained from molecular 





dynamic simulations using microscopic principles. 


B. BACKGROUND 








Predicting the physical properties of macroscopic 
systems from the microscopic nature of their constituent 


atom is a difficult problem. Yet the problem is one that 





must be addressed if we are to develop useful technologies 


based on embedded nano-scale sensors. Traditionally there 





have been two ways of tackling this kind of problem. The 


first is a first principles, “brute force” method using 





quantum mechanics. This method is severely constrained by 
the size of the system it can handle. A truly first- 


principles approach in predicting the behavior of systems 





where quantum effects are significant, is prohibitively 








costly in terms of computing power. The second method is 
molecular dynamics, which solves the problem by making 
specific assumptions from the onset to simplify the 
problem[3]. 


In this thesis we will see that a combination of 





techniques borrowed from molecular dynamics and quantum 
mechanics can be used to model nano-structures. Both 
molecular dynamics and quantum mechanics will be discussed 


separately below, followed by a discussion on how they 





merge. 


1. Molecular Dynamics 





In almost all cases, molecular dynamics simulations 





are based on the Lennard-Jones potential of inter-atomic 


interactions. For a pair of atoms i and j located at Ff and 


r the potential: énergy 1s given ‘by [41% 


where F,=Fr-Tr,, r=|7| and € and tT are phenemenological 


parameters specific to each material system. The Lennard- 


Jones potential features a short-range repulsive potential 





together with a long-range attractive potential. The 1/r‘° 
attractive potential arises from fluctuating dipole moments 


on neutral atoms, while the l/r’ term models the short-range 


repulsion that arises from the Pauli exclusion effects of 











overlapping electron clouds. The Lennard-Jones potential 


defines the repulsion, attraction and finally the cutoff 








range where the potential no longer is applicable. The 








first simplification in molecular dynamics is to ignore the 








attractive tail. The above equation can be modified to: 
12 6 
- T T _ 41/6 
u(r, )=4¢ —} -|— toe Tes = 2 [ 
r r 
ij ij 
with r, chosen such that u(r.)=0. What we want is the force 





of the interaction that is given by f =-Vulr), which results 


in: 


= 48 ]( +t oy ines I 


re 
y 2 y 
T r.. 2\ 7, 





As r aneréeases to +, the force drops to zero, which: prevents 


any discontinuities. 


So far only a two-atom system has been discussed. For 








an N-atom system, the total intermolecular potential 


3 





FUnCELOn 1s: represented by u(r”). Relsbing bhi Oo: Our 


Lennard-Jones pair potential u(r) we have foi 


ur”)=V YuG,). 


Jj ois<j 





Since, at the micro-scale, there are no dissipative forces 





acting among molecules, intermolecular forces are 


conservative and the total force acting on molecule i is: 


which is equivalent to: 








Therefore the kinetic and potential energies, per 


atom, can be written as: 


4C Tin 8 6 
E hot = oe Mees ae ae 
oe N pis ij ij 
where E,,=E,,,+£,, as usual. 


kin pot 


Molecular dynamics simulations compute the motions of 





individual molecules. From this motion, molecular dynamics 





Simulation has direct access to the energy and pressure of 


a system in equilibrium. Energy and pressure can be 





expressed in terms of temperature and density. It is the 
measurements of the temperature and density of a system, 
using molecular dynamics, which provides access to the 


thermodynamics of a system and it is the thermodynamics of 





a system that give us a detailed description of the 


microscopic behavior. 





If we turn this around and specify the position of the 





atoms of a system, the energy can be determined. This 


4 








scheme is called th mbedded atom method (EAM)[6]. Daw 
[7] also showed that the energy of a system of atoms is 
uniquely specified by the its electron density. From the 


system energy, by using numerical methods, the forces can 





be determined as well as the thermodynamic properties of 


the material. The embedded atom method is one schemes used 





in order to obtain system energies. The energy of a system 





can be represented again by [8]: 


Fess = Ekin as Eso , OF 


1 
Ee = LF a (p,,; )+ eos (r,,) ' 


where p,; ts the total electron density seen by atom i due 





to the rest of the atoms in the system, E, is the embedding 


ee 








energy for placing an atom into that electron density and 





0, is the short range pair interaction representing the 





core-core repulsion of atoms i and j separated by distance 
rj- Quantum mechanic provides a means of obtaining these 
quantities. 

2: Quantum Mechanics 





In this section, only highlights of the important 


parts will be discussed. 





In a perfect world, the exact solution for the time 


independent Schrodinger equation for an N-electron wave 





function could be determined. However determining the 





exact solution is complex if not “hopeless”[9]. Basically 





the problem is to solve the Schrodinger wave equation in 


its time dependant form for a particle of energy moving in 





a potential V in three dimensions[10]: 
h?>(vov ov ow 
a a Ye ae 
2m|\| ox oy dz 














Jreveonede nt, 
ot 


where W(x, y,z) is the wave function for N atoms and 


potential Vcould be a function of space and time. Using 





variational (Hartree-Fock method) calculations of self- 
consistent fields, one can narrow down a solution. By 


incorporating the Slater determinant, one can identify a 











solution that includes the requirements of the Pauli 








exclusion pricipal. A direct, frontal assault on the many 





electron Shrodinger equation is impractical for systems 
containing more than a few electrons. We are guided, 
however, in developing approximate methods by the Ritz 


variational principle, which says[1ll]: 


Waly) 
oe Sy) 





where H is the Hamiltonian and is given by[12]: 


N 2 2 2 
‘ 2 7, 
H=) pent caer Wo 


: rd fale 





By taking into account the normalization[13] 

Jarig@)P=1 

and with some manipulation, we end up with the Hartree-Foch 
equations. The problem is that the Hartree-Fock method 
works for atoms with few electrons. Even with today’s 


computers this method is cumbersome. 


The Thomas-Fermi method, introduced in the 1920's, 





showed promise in evaluating the energy for atoms with many 





electrons. In this method the potential is expressed in 


terms of the particle number density at a position vice the 





traditional probability density, as used above. For 


example[14]: 


1 2 
n(r)= as 2mV(r))2 . 


The above equation is a semi-classical approximation that 





assumes that there are many electrons within volume 
elements in which the potential is nearly constant and the 


wavelength is small since the electrons are assumed to be 





in a state of high energy[15]. The problem with the 








Thomas-Fermi method is that it is not very accurate for 











small ranges, since the potential changes become large. 








Also it is not very accurate for large ranges when th 


wavelength becomes large. 





Density Functional Theory changes the form of the 


solution for which we are looking. Specifically, given a 








potential, the above solution results in an electron 








density, where th nergy is a functional of that 








potential. Hohenberg and Kohn[16] stated that there 





existed a functional defined for all non-degenerate ground 





state densities such that, for a given potential, the 


quantity has its unique minimum for the correct ground 








state density, associated with the potential. In other 


words, the energy is a function of the density, and the 





potential can be determined by the electron density (to 





within a constant). 


Next Stott and Zaremba[17][18] proved a useful 





corollary, in that the energy of an impurity in a host is a 
function of the function (or functional) of the electron 


density of that same host without impurities, or 


E, =F [p,(r)l, where p, (r) Ls the pure host electron density. 
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From the corollary, the conclusion was drawn that the 


embedding energy of an impurity can be determined from the 





electron density of the host before the impurity is added. 








Daw and Baskes[19] then showed that the energy then is a 














function of the electron density at the site of the 








impurity plus an additional electrostatic interaction term. 





Mathematically: 
1 
1 Oe S F,(p, ()+>Lol), 
j 


where: (Pp; Us the electron density of the host. without 
impurity, and @® is the short range electrostatic pair 


potential. The total energy of a system is the 


superposition of all contributions[20]: 


Ew =D FulPuln))+> L(y) 





3. Paradyn 
Paradyn[21] is the computer code used in this thesis 


to solve numerically for the total energy of a metallic 





based system. From the code, thermodynamic properties are 





obtained that can be used for a variety of purposes. The 


properties required by this thesis are the global off 





diagonal stresses and the elongation length of a cubic 





metallic lattice, due to a tensile loading stress. 





There are two main input files required to run 





Paradyn. The first file contains the embedded potential, 





and pair interaction information. The second file sets up 





the thermodynamics for determining the kinetic energy. 
The potential file contains the following information: 
atomic number; atomic mass, lattice spacing, potential 


cutoff range for each atom, grid of inter-atomic force, 


grid of atom electron density and the spherically averaged 
atomic density[22]. The grid of the atom electron density 
is interpolated for each atom site in order to determine 


the electron density at each atom. The inter-atomic force 





is also interpolated for each atom site and multiplied with 








the atom density in order to determine th mbedding energy 





of each atom. The embedding energy is then summed over all 





the host atoms keeping to the predetermined cutoff range, 


in order to determine the cohesion energy of the system. 








Finally, the electron density at each atom is multiplied by 














th mbedded energy, which determines the force on each 


atom. The components of the individual atomic force 





tensor, multiplied by the inverse of the corresponding 








neighbor range squared, are used to determine the stress 
tensor for each atom. 


The thermodynamic file takes the temperature required 








by the user and computes the local velocity of each atom by 








v=~vV2Tk/m. In order to ensure the system represents a real 





model, random numbers are generated and assigned to the 


initial velocities, which are used to develop the overall 

















velocity distribution, so that v,=v*randomnumber.  v, and v, 
are determined in the same fashion. The kinetic energy of 
each atom is determined from the thr velocity components 

















m 
Bag Wet, ete) The individual kinetic energies are then 

















summed over all the atoms to determine the total kinetic 
energy of the system. The user using this file also inputs 
the number of atoms that define a system. The velocities 


are also used in determining the diagonal components of the 


stress tensor for each atom. For example: 6, =mv,v,). Once 








the stress tensor is obtained the pressure of the system 


can be determined by: 





P= ( prot, )/3Vol 


With the stress tensor developed by both the kinetic 





energy and potential energy routines, the velocity of the 





kinetic energy routine is adjusted in an iterative loop 





until the temperatures of the two separate routines are 











within tolerance. Depending on the initial user flags, the 





tolerance can be dictated by the user. For the program 





details, refer to either the Paradyn computer code or 








supporting documentation. 


4. Elastic Properties 





In general the microscopic stress tensor is of the 


form: 
a xy J. 
J=|J re J + J a 
J, J 








where Jj measures the rate at which the B-directed 
momentum is transported in the Q@-direction[23]. The stress 
tensor components are made up of both kinetic and potential 
energy parts as shown below: 


N 1 N 
Jap = mY Vialip eee Tig Fie. 


i#j 





where the first term is equivalent to the kinetic energy 





and the second term the potential energy. For a homogenous 





isotropic material the stress tensor is symmetric. In 








order to obtain the stress auto-correlation function, the 
off diagonals are used in the form: 


10 


nl) = Le es to) Zap (fo +9) 


Joep 2S the ‘sum “of two. terms, thus: the: corre lation 
function N(t) has three contributions: (a) a 
kinetic term, which measures the correlation of 
the momentum transport caused by atomic motions; 











(b) a potential term, which measure the 
correlation of momentum transport cause by inter- 
atomic forces; and (c) a cross term which 





measures the coupling of atomic motions and 
forces. [24] 





For programming purposes the above equation can be 


simplified to: 


n)= a NM Lbs op (ty ap (tn +1) 


where M is the number of samples. In the static limit as 








t goes to 0, the shear modulus n(0) can be obtained. This 





thesis compares the shear modulus obtained by microscopic 





means and compared to that obtained for macroscopic systems 


to ensure the modulus is preserved in some form or another. 
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II) RESULTS AND DISCUSSION 


This section will outline the procedure used in the 





determination of the Modulus of Elasticity and the Shear 


Modulus from the data produced from multiple Paradyn 








Simulations. The mechanical properties will then be 





directly compared to known values. Five substrate 


materials were selected for this baseline comparison. 











Specifically copper, nickel, gold, silver and platinum and 





monel. It should be noted that all references to the auto- 
correlation function refer to the C++ program at Appendix 


A. 


A. MODULUS OF ELASTICITY 





Initially Paradyn was run with no tensile load. The 








no load case provided a baseline of how the lattice 


structure equilibrated and subsequently vibrated. The 








vibrations are due to the kinetic energy of the individual 


atomic structures. Once the no load case data was 





obtained, the sample was subjected to a load of 10 MPa of 


tension. Again the sample increased in size and then 





oscillated. For this thesis the structure was allowed to 


increase to its equilibrated size, then oscillate through 








two full periods. In Figure 1, a sample of platinum was 








allowed to equilibriate at loads of 0 MPa and 10 MPa of 


tension. The initial expansion of the sample is evident. 





Once the initial expansion is over the sample can be seen 








to oscillate. In Figure 2, the oscillations over two full 


periods are expanded. The increase in separation of the 


1.3 








equilibrium length is what is critical in determining the 





elastic modulus. 


Longitudinal Oscillations in a Platinum Sample 


= 
oO 
a 


= 
oO 
i) 


——-0 MPa Tension 
1 MPa Tension 


-~ 
<x 
— 
= 
i 
D 
ra 
® 
_ 
£ 
® 
D 
ro 
© 
= 
oO 


400 800 


Number of Time Steps 








Figure 1. Longitudinal Changes in Length of a Platinum 
Sample With No Load and 10MPa of Tension. 
Temperature is 300K. One Time Step is 5 ps. 





Longitudinal Oscillations in a Platinum Sample 
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Change in Length (A) 
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Figure 2. Expansion of Figure 1. The Longitudinal 
Changes in Length of a Platinum Sample With No 
Load and 10 MPa of Tension. Expansion Consists 





of the Re-centered Time Steps 110-540 From Figure 
1, In Order to Encompass Two Full Periods of 
Oscillations. 
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In determining the elastic modulus, first the non- 
loaded sample was allowed to expand and come to equilibrium 
for a minimum of two complete oscillatory periods. The 
average amplitude of the period peaks was used to determine 


the mean length of the sample (J,). Next a tensile stress 





(O) was applied to the sample. Again the sample was 
allowed to expand and cycle through at least two full 


cycles. The mean length of the period peaks was used to 








determine th longated length l. The strain (€) was then 





determined by e€=(1,—-J)/1,. The modulus of elasticity (EF) 
was then calculated in the normal way by E=o/€. The 
modulus of elasticity determined from Paradyn is listed and 


compared in Table 1. 





Modulus of Elasticity (GPa) lus of Elasticity (GPa) 
Theoretical from Reference (25) Calculated from Paradyn 


220.7 252.0 
146. 





161.7 








4 3 


Table 1. Comparison of the Theoretical and 
Calculated Modulus of Elasticity 





A sample of the alloy monel was also tested. The 





sample of monel consisted of 64% nickel atoms and 33% 
copper atoms, where the copper atoms were placed ina 
nickel substrate at random. The result for monel is as 
expected. The modulus of elasticity was expected to be 


between the values obtained for pure samples of nickel and 




















copper. From Table 1, the simulated modulii of elasticity 
are greater than the theoretical values. It is the 
difference in methodology between the techniques that plays 
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an important role. The modulli resulting from the 


simulations are based on pure metallic samples, whereas the 





theoretical values are based on bulk material test samples. 








These test samples are made of real materials that contain 





impurities and crystalline structures. 


The copper sample was exposed to increasing 


temperatures and increasing tensile stresses. The results 





of the increasing temperature runs are found in Table 2 and 


the results of the increasing tensile stress are found in 





Figure 3. 





Average Internal Temperature Modulus of Elasticity 


(K) (GPa) 


220.8 


220.8 


Table 2. Modulus of Elasticity as a Function of 
Temperature for a Sample of Pure Copper. Results 
are Based on Simulation and are within 15% of the 
theoretical value. The Melting Point of Copper 

is 1084.62 C [26] 
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Figure 3. Stress Strain Curve of a Sample of Pure 
Copper. Data Obtained from Simulation. 
Theoretical Yield Strength of Copper is 27.5 MPa. 
[27] 


Further analysis of the Elastic Modulus was completed 


using the same copper sample, in that dislocations were 





added to the sample. In the first run 10 copper atoms were 
removed at random, run 2 had 20 atoms randomly removed and 
finally a block of 28 atoms were removed from the center of 


the sample. The results of these runs are found in Table 3 


Description Applied Tensile odulus of 
—— Elasticity 
(GPa) 











500 Atom Sample 


_ Random Atoms Removed 

Random Atoms Removed 

= _s (Block of 28 Atoms) Removed 

Table 3. Modulus of Elasticity for a Sample of 
Pure Copper with Dislocated Atoms. Temperature 





























Constant at 300K, With an Applied Tensile Stress 
of 10 Bar. 


For all the runs in which the Modulus of Elasticity 


was determined, the values were comparable to the 
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theoretical values. From at Table 1, where the six 








baseline substrates are directly compared, the modulus of 
elasticity for these six materials are not only in the 


correct order but they have values within 15% of their 








respective theoretical values. 





Table 2 illustrates that as the temperature of the 
copper sample is increased, the modulus of elasticity 


decreases. This trend is correct and is one reason that 





many metal materials are worked at higher temperatures. It 





takes less energy to form the material at the higher 








temperature, because it is less stiff. 


Figure 1 is the result of a plot of stress versus 
strain for the copper sample. The plot of the resulting 


Paradyn output data is representative of a typical Stress-— 








Strain curve plot of a metallic material. 
For the copper sample with dislocations, Table 3, 
illustrates that as atoms are removed th lastic modulus 








decreased as expected. Also, there is a difference in the 


modulus of elasticity in how the atoms were removed. For 








random atom removals, the modulus decreases gradually as 
more and more atoms are removed. For the removal of a 


block of atoms, the modulus decreases abruptly. Of note, 





the run data show that the atoms migrate towards the 


dislocations and take up new equilibrium positions as 





compared to the sample with no dislocations. As well, for 





the runs where 20 or more atoms were removed, the sample 





volume decreased overall. Although this decrease in volume 


is in the order of cubic angstroms, it was none the less 





present. 
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B. SHEAR MODULUS 





As described earlier, molecular dynamic simulations 
typically compute the shear modulus by the virial 
expression: 


(ata St ue) 


Im n\t)= 
10 3kT NM = 








A Ct+ program was written in order to determine th 


shear modulus from the Paradyn stress output file. The 





code for this program is included at annex A. 


The theoretical value of the shear modulus was 





determined by[28]: 


E 
G = ——; 
2(1+y) 
where E is the modulus of Elasticity and y is the 


material's Poisson ratio. 


Table 4 contains the calculated theoretical shear 


modulus and the shear modulus computed from the stress 





auto-correlation function based on Paradyn output data. 














Poisson Ratio Theoretical Shear Experimental Shear 
[25]*, [2917[30]* Modulus (GPa) Modulus (GPa) 
pe $5258 — +] = 54.5 ———. 
aa i eee ee ro eee 
Table 4. The Theoretical Shear Modulus as 
Determined from the Modulus of Elasticity (given 
in Table 1) and the Poisson Ratio. The 
theoretical and experimental shear modulus' are 





presented side by side for direct comparison. 


From Table 4, the resulting shear modulii computed 


from Paradyn data are representative of the theoretical 
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values. For nickel, platinum and copper, the values of the 





shear modulus agree to within 15%. For silver and gold the 


experimental shear modulus is off by approximately 30-50% 
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III CONCLUSIONS AND RECOMMENDATIONS 


A. CONCLUSIONS 


Bulk properties of materials in macroscopic systems 





are well defined and can easily be obtained from 





destructive testing of material samples. Computer 





Simulations have also been successful at modeling these 





properties on a macroscopic scale. As technology advances 





to produce sensors and machines on a nano-scale, the task 





at hand is to determine the mechanical properties of these 


micro-systems so that when used they do not loose these 





properties and become a single point of failure ina 


critical system. Therefore it is important that a good 








model would combine Molecular Dynamics and Quantum 








Mechanics in order to determine the mechanical properties 





from interactions on an atomic level. 








In this thesis the Paradyn model was used to determine 
the mechanical properties of pure metallic samples 
containing about 500 atoms. Preliminary analysis of the 


data from the Paradyn simulations is encouraging. The 





results show that, for the six substrate samples chosen, 





the modulus of elasticity can be obtained to within 15% of 

















the theoretical value. In the simulations where the copper 
sample was exposed to increasing tensile stress, 
temperature and number of dislocations, the trend for the 











elastic modulus is reasonable. In plotting the stress and 





strain for a copper sample, the resulting stress strain 


curve was shown to be a typical representation of a 





metallic material. Finally, the simulation where atoms 
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were dislocated from the sample material illustrates that 
Paradyn was capable of producing results that were 


reasonable. 


From the resulting data it can be concluded that 





Paradyn has the ability to simulate a micro-metallic system 
from which the bulk mechanical properties of the system are 


preserved. These results also indicate that the Paradyn 











model does accurately determine the inter-atomic forces 


between atoms of a metallic sample. 


B. RECOMMENDATIONS 





This thesis provides baseline material for obtaining 
mechanical properties of pure microscopic metallic systems 


using the Paradyn model of molecular dynamic simulation. 








In order to advance research into incorporating the 
embedding of micro technology to include carbon nanotube 


structures the following recommendations are made. 





Conduct simulations for alloy and composite 
materials that could also be used as potential 


substrate materials; 





2% Determine if the model accurately reflects 





material defects such as hydrogen embrittlement 


or atomic impurities. 





3s Further research the determination of the shear 





modulus in order to determine if the model does 





accurately compute the microscopic stress tensor 





from the inter-atomic forces. 


a. Conduct simulations where carbon nanotube 


structures ar mbedded into substrates. 





22 


APPENDIX A C++ CODE TO DETERMINE THE SHEAR MODULUS 


[ [KR RRR KK RK KR KR KK KR KK eK eK KK KK IK Kk 


//File: stress.cpp 

//Name: Richard Duff 

//PHOO80 Thesis 

[fi March. 2003 

//Description: This program calculates the modulus 
// by using the off diagonal relaxation stesses 


// from paradyn and solving the using autocorrelation 
[ [KR KKK RK RK KR KR IK eK eK RK eK eK KK KK IK Hk 











#include <stdio.h> 
#include <math.h> 
#include <stdlib.h> 














// Initial values 








double boltz=8.62E—-5, msteps=20.0, natoms=500.0; 
double temp=147.5784, vol=5452.596; 











//Variables 
// 
// 
AE. As gu See 

float fa[25], fb[25], fel[25]; 
double ACF, etat, rho; 





main () 


{ 





//Initialize setup parameters 
rho=natoms/vol; 

//Set up for file read and write 
FILE *fp; 





//read stress matrix 
fp=fopen("ninonpotstress.txt","r"); 
for (1=0;1<20;++1) 

{ 





fscanf(fp, "SE SE SE", &fali]l, &fb[i], &fcli]) 


} 
fclose(fp); 
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// this. loop 18 only for verification pirposés to enstre 
input is correct 

fp=fopen("ninonpotcopy.txt","w"); 
fprintf(fp,"Input - The Off Diagonal Stress Tensor 
Components Are:\n"); 

EPEINCE (tp 7% WilZ J13 
J23\n\n"); 

for (1=0;1<20;++4+1) 

{ 

fprintf (fp,"%14.5E $14.5E $14.5E\n", fa[i], fb[il, 
BO La) 7 














} 
fclose(fp); 
// Open file for output data 
fp=fopen ("ninonpotmodulus.txt","w"); 
PrintEe(™t ACE) % 





//loop for number of time origins 
fprintf(fp,"Nickel room temp not minimized using pot 
file only"); 





fprintf(fp,"\n t ACF 
Eta(t)\n\n"); 
for (t=0;t<20;++4+t) 


{ 
// loop for number of sample points 
ACF=0.0; 
for (i=0;1i<(20-t);++i) 
af 


ACF=ACF+ (fa[i]*fa[itt]+fb[i]*fo[itt]+fcli] *fclitt]); 
} 
ACF=ACF/ (20); 
printf ("%S4d $14.5E\n",t,ACF); 
etat=(rho*ACF) /(3*boltz*temp*natoms) ; 
fprintf(fp," %4d %14.5E 
$14.5E\n",t,ACF,etat); 





} /fend t 
fclose(fp); 
return (0); 


} //end main () 


//end of file stress.cpp 
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